####################################################################
# Part 4: SIR Algorithm to Approximately Sample from the Posterior #
####################################################################
# April 30, 2024

# Npostsamp = size of sample generated by the SIR algorithm
Npostsamp = 10

# This function generates an approximate iid sample of Npostsamp from the posterior using the SIR algorithm

SIR_algorithm = function(Npostsamp, cum_weights, p, mu_xi, xi, Sigma){
  #' Using the SIR algorithm from Rubin, values of i in {1, 2, ..., Npostimp}  
  #' and corresponding values of mu, xi and Sigma.   
  #' @param Npostsamp the size of the sample that will be generated.
  #' @param cum_weights the vector containing the cumulative weights.
  #' @param mu the mu matrix.
  #' @param xi the list containing xi's.
  U = runif(Npostsamp)
  i = findInterval(U, cum_weights)+1
  sample_mu_xi = mu_xi[i, ]
  sample_xi = xi[,,i]
  sample_Sigma = Sigma[i]
  newlist = list("sample_mu_xi" = sample_mu_xi, "sample_xi" = sample_xi, "sample_Sigma" = sample_Sigma)
  return(newlist)
}

post_SIR_sample = SIR_algorithm(Npostsamp, cum_weights = imp_cdf, 
                                p = p, mu = imp_mu, xi = imp_xi, Sigma = imp_Sigma)


mu_postsamp = post_SIR_sample$sample_mu_xi
xi_postsamp = post_SIR_sample$sample_xi
Sigma_postsamp = post_SIR_sample$sample_Sigma